#!pip install skforecast
#Librerías para cálculos y gráficos
#====================================================================================
import pandas as pd
import numpy as np
from numpy import array
import plotly.graph_objs as go
import matplotlib.pyplot as plt
#pip install matplotlib==3.7.1
import seaborn as sns
#Librerías para Preprocesamiento de datos
#====================================================================================
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import OneHotEncoder
from sklearn.compose import ColumnTransformer
#Librerías para Implementación de Modelos
#====================================================================================
from sklearn.ensemble import RandomForestRegressor
from skforecast.ForecasterAutoreg import ForecasterAutoreg
from skforecast.ForecasterAutoregCustom import ForecasterAutoregCustom
from skforecast.model_selection import grid_search_forecaster
from skforecast.model_selection import backtesting_forecaster
from sklearn.metrics import mean_squared_error, r2_score, mean_absolute_error
from skforecast.utils import save_forecaster
from skforecast.utils import load_forecaster
Importamos la base de datos con la cual vamos a trabajar
dfmacro = pd.read_csv(r'C:\Users\Vivi\Downloads\BaseAccionesF.csv')
dfmacro.columns = ["Date","GSPC","AAPL","AMZN","MSTF","TSLA","GOOG","GOOGL","NVDA","BRK.B","META","UNH","JNJ","PG","VIX","UNRATE","FEDFUNDS","UMCSENT","CPIAUCSL","DOLARINDEX"]
dftotal = pd.read_csv(r'C:\Users\Vivi\Downloads\BaseAcciones2.csv')
dftotal['Date'] = dfmacro['Date']
dftotal = dftotal.dropna()
dftotal.iloc[:,[43,56,75]].describe()
| GSPC.ReturnSuavizado | VIX.ReturnSuavizado | DolarIndex.ReturnSuavizado | |
|---|---|---|---|
| count | 2667.000000 | 2667.000000 | 2667.000000 |
| mean | 0.062439 | -0.344104 | 0.008501 |
| std | 0.880420 | 7.081882 | 0.298783 |
| min | -3.002265 | -26.622753 | -1.188374 |
| 25% | -0.360074 | -4.487971 | -0.166130 |
| 50% | 0.063567 | -0.749383 | -0.001720 |
| 75% | 0.544763 | 3.449533 | 0.181680 |
| max | 3.101521 | 26.760500 | 1.222749 |
df2= dftotal.loc[:,['Date','GSPC.ReturnSuavizado', 'AAPL.ReturnSuavizado', 'AMZN.ReturnSuavizado',
'MSTF.ReturnSuavizado', 'TSLA.ReturnSuavizado', 'GOOG.ReturnSuavizado',
'GOOGL.ReturnSuavizado', 'NVDA.ReturnSuavizado',
'BRK.B.ReturnSuavizado', 'META.ReturnSuavizado', 'UNH.ReturnSuavizado',
'JNJ.ReturnSuavizado', 'PG.ReturnSuavizado', 'VIX.ReturnSuavizado','DolarIndex.ReturnSuavizado']]
df2['Date'] = pd.to_datetime(df2.Date)
df2 = df2.set_index('Date')
df2 = df2.asfreq('B',method='bfill')
#df2 = df2.iloc[:,1:]
print(f'Número de filas con missing values: {df2.isnull().any(axis=1).mean()}')
Número de filas con missing values: 0.0
# Verificar que un índice temporal está completo
# ==============================================================================
(df2.index == pd.date_range(
start = df2.index.min(),
end = df2.index.max(),
freq = df2.index.freq)
).all()
True
# Verificar que no existan valores nulos
# ==============================================================================
print(f'Número de filas con missing values: {df2.isnull().any(axis=1).mean()}')
Número de filas con missing values: 0.0
# Definir tamaños de train y test
# ==============================================================================
dataindex2= df2.index
steps = 20
n100 = len(df2)
ntrainall = int(n100-steps)
ntestall = n100 - ntrainall
ntrain = int(ntrainall-steps)
ntest = steps
n100, ntrainall, ntestall,ntrain,ntest
(2812, 2792, 20, 2772, 20)
# Construir conjuntos de entrenamiento y validación
# ==============================================================================
dftrainall = df2[:ntrainall]
dftestall = df2[ntrainall:]
len(dftrainall),len(dftestall)
(2792, 20)
# Sacamos una porción de datos de train para validación nuevamente
# ==============================================================================
dftrain_sinest = dftrainall[:ntrain]
dftest_sinest = dftrainall[ntrain:]
dftestall_sinest = df2[ntrainall:]
len(dftrain_sinest),len(dftest_sinest),len(dftestall_sinest)
(2772, 20, 20)
scaler = StandardScaler()
cols_to_standardize = list(dftrain_sinest.columns)
dftrain = pd.DataFrame(scaler.fit_transform(dftrain_sinest[cols_to_standardize]),
index=dftrain_sinest.index, columns=cols_to_standardize)
# Estandarizamos los datos de test
dftest = pd.DataFrame(scaler.transform(dftest_sinest[cols_to_standardize]),
index=dftest_sinest.index, columns=cols_to_standardize)
dftestall = pd.DataFrame(scaler.transform(dftestall_sinest[cols_to_standardize]),
index=dftestall_sinest.index, columns=cols_to_standardize)
#dftrain.to_excel('dftrain.xlsx')
#dftest.to_excel('dftest.xlsx')
#dftestall.to_excel('dftestall.xlsx')
En primera instancia se corre el modelo con el 100% de las variables con las que contamos, dejando a y como predictor y en variables exógenas los rendimientos de las acciones y las variables macroeconómicas.
Se calculan los hiperparámetros y se le pide al modelo que retorne el modelo mejor entrenado para las predicciones:
# Cálculo de hiperparámetros por grid search
# ==============================================================================
forecaster1 = ForecasterAutoreg(
regressor = RandomForestRegressor(random_state=123),
lags = 14
)
# Lags used as predictors
lags_grid = [5, 20]
# Regressor's hyperparameters
param_grid = {'n_estimators': [100, 500],
'max_depth': [3, 5, 10]}
results_grid = grid_search_forecaster(
forecaster = forecaster1,
y = dftrain.iloc[:,0],
exog = dftrain.iloc[:,1:],
param_grid = param_grid,
lags_grid = lags_grid,
steps = steps,
refit = True,
metric = 'mean_squared_error',
initial_train_size = int(len(dftrain)*0.5),
fixed_train_size = False,
return_best = True,
verbose = False
)
Number of models compared: 12.
loop lags_grid: 0%| | 0/2 [00:00<?, ?it/s] loop param_grid: 0%| | 0/6 [00:00<?, ?it/s]C:\Users\Vivi\tf\lib\site-packages\skforecast\model_selection\model_selection.py:370: RuntimeWarning: The forecaster will be fit 70 times. This can take substantial amounts of time. If not feasible, try with `refit = False`. warnings.warn( loop param_grid: 17%|██████▎ | 1/6 [01:27<07:18, 87.69s/it] loop param_grid: 33%|████████████▎ | 2/6 [08:36<19:13, 288.40s/it] loop param_grid: 50%|██████████████████▌ | 3/6 [10:45<10:46, 215.57s/it] loop param_grid: 67%|████████████████████████▋ | 4/6 [21:24<12:45, 382.61s/it] loop param_grid: 83%|██████████████████████████████▊ | 5/6 [25:06<05:24, 324.87s/it] loop param_grid: 100%|█████████████████████████████████████| 6/6 [43:20<00:00, 586.36s/it] loop lags_grid: 50%|██████████████████▌ | 1/2 [43:20<43:20, 2600.63s/it] loop param_grid: 0%| | 0/6 [00:00<?, ?it/s] loop param_grid: 17%|██████▏ | 1/6 [02:18<11:34, 138.91s/it] loop param_grid: 33%|████████████▎ | 2/6 [13:50<30:55, 463.84s/it] loop param_grid: 50%|██████████████████▌ | 3/6 [17:24<17:30, 350.02s/it] loop param_grid: 67%|████████████████████████▋ | 4/6 [35:14<21:08, 634.33s/it] loop param_grid: 83%|██████████████████████████████▊ | 5/6 [41:31<09:01, 541.35s/it] loop param_grid: 100%|███████████████████████████████████| 6/6 [1:12:48<00:00, 995.51s/it] loop lags_grid: 100%|███████████████████████████████████| 2/2 [1:56:09<00:00, 3484.60s/it]
`Forecaster` refitted using the best-found lags and parameters, and the whole data set:
Lags: [1 2 3 4 5]
Parameters: {'max_depth': 10, 'n_estimators': 500}
Backtesting metric: 0.3251563573906958
# Predicciones
# ==============================================================================
predicciones = forecaster1.predict_interval(steps=steps, exog=dftest.iloc[:,1:],interval = [1, 99],n_boot = 500)
#Podemos primeramente graficar las secciones: la serie de entrenamiento, la predicción en la sección de prueba y la sección de prueba.
def grafico_comparativo(dftrain,dftest,predicciones):
trace1 = go.Scatter(
x = dftrain.index,
y = dftrain.iloc[:,0],
mode = 'lines',
name = 'Data'
)
trace2 = go.Scatter(
x = dftest.index,
y = dftest.iloc[:,0],
mode = 'lines',
name = 'Prediction'
)
trace3 = go.Scatter(
x = predicciones.index,
y = predicciones['pred'],
mode='lines',
name = 'Reality'
)
layout = go.Layout(
title = "S&P Pronóstico",
xaxis = {'title' : "Date"},
yaxis = {'title' : "Close"}
)
fig = go.Figure(data=[trace1, trace2, trace3], layout=layout)
fig.show()
grafico_comparativo(dftrain,dftest,predicciones)
def grafico_intervalos(y,predicciones):
fig, ax = plt.subplots(figsize=(10, 3))
y.plot(ax=ax, label='test')
predicciones['pred'].plot(ax=ax, label='predicciones')
ax.fill_between(
predicciones.index,
predicciones['lower_bound'],
predicciones['upper_bound'],
color = 'red',
alpha = 0.2
)
ax.legend();
grafico_intervalos(dftest.iloc[:,0],predicciones)
def grafico(y,pred_1fin):
fig, ax = plt.subplots(figsize=(10, 3))
y.plot(ax=ax, label='test')
pred_1fin['pred'].plot(ax=ax, label='predicciones')
ax.legend();
def desempeño(dftest,predicciones):
error_mse = mean_squared_error(
y_true = dftest.iloc[:,0],
y_pred = predicciones['pred']
)
rmse = np.sqrt (error_mse)
print(f"EL MSE del modelo en test es de: {error_mse}")
print(f"EL RMSE del modelo en test es de: {rmse}")
r2 = r2_score(
y_true = dftest.iloc[:,0],
y_pred = predicciones['pred']
)
print(f"EL R2 del modelo en test es de: {r2}")
# suponga que y_true son los valores reales y y_pred son las predicciones del modelo
error_mae = mean_absolute_error(
y_true = dftest.iloc[:,0],
y_pred = predicciones['pred']
)
print(f"EL MAE del modelo en test es de: {error_mae}")
desempeño(dftest,predicciones,)
EL MSE del modelo en test es de: 0.13204347393516816 EL RMSE del modelo en test es de: 0.36337786660055144 EL R2 del modelo en test es de: 0.9011377899512703 EL MAE del modelo en test es de: 0.36337786660055144
Con los datos obtenidos se corre el modelo sugerido con el 100% de los datos de entrenamiento y se obtiene:
dftrainall_est = pd.concat([dftrain,dftest])
# Configuramos la función de pronóstico de forecaster con la mejor combinación obteida de hiperparámetros.
# ==============================================================================
### Aquí debe ir tu codigo,
regressor = RandomForestRegressor(max_depth= 10, n_estimators=500 , random_state=123)
### Hasta aquí modificas
forecaster = ForecasterAutoreg(
regressor = regressor,
lags = 5
)
forecaster.fit(y = dftrainall_est.iloc[:,0], exog = dftrainall_est.iloc[:,1:])
# Pronóstico, en este caso elegimos una ventana de 20 días.
# ==============================================================================
pred_1fin = forecaster.predict_interval(steps=steps, exog=dftestall.iloc[:,1:],interval = [1, 99],n_boot = 500)
grafico_intervalos(dftestall.iloc[:,0],pred_1fin)
desempeño(dftestall,pred_1fin)
EL MSE del modelo en test es de: 0.12163476353394123 EL RMSE del modelo en test es de: 0.3487617575565607 EL R2 del modelo en test es de: 0.9082587299062685 EL MAE del modelo en test es de: 0.25495078018323786
forecaster.get_feature_importance()
| feature | importance | |
|---|---|---|
| 0 | lag_1 | 0.008050 |
| 1 | lag_2 | 0.005388 |
| 2 | lag_3 | 0.005674 |
| 3 | lag_4 | 0.006963 |
| 4 | lag_5 | 0.006396 |
| 5 | AAPL.ReturnSuavizado | 0.043067 |
| 6 | AMZN.ReturnSuavizado | 0.023042 |
| 7 | MSTF.ReturnSuavizado | 0.176339 |
| 8 | TSLA.ReturnSuavizado | 0.009582 |
| 9 | GOOG.ReturnSuavizado | 0.047308 |
| 10 | GOOGL.ReturnSuavizado | 0.044740 |
| 11 | NVDA.ReturnSuavizado | 0.026820 |
| 12 | BRK.B.ReturnSuavizado | 0.328805 |
| 13 | META.ReturnSuavizado | 0.008780 |
| 14 | UNH.ReturnSuavizado | 0.016280 |
| 15 | JNJ.ReturnSuavizado | 0.013042 |
| 16 | PG.ReturnSuavizado | 0.014296 |
| 17 | VIX.ReturnSuavizado | 0.204769 |
| 18 | DolarIndex.ReturnSuavizado | 0.010659 |
Al analizar la importancia de los predictores, encontramos que el VIX, MSTF y BRK.B son los 3 predictores de mayor importancia. Dado que las variables que ingresamos como exógenas al modelo requieren ser conocidas y que esto implicaría tener un modelo predictivo por variable, vamos a explorar la forma de reducir la cantidad de variables que tenemos que predecir.
Analicemos estas dos acciones dentro del portafolio global:
import seaborn as sns
fig, ax = plt.subplots(figsize=(10,6))
n = len(df2.columns)
df_cor = df2.reset_index()
df_cor = df_cor.iloc[:,1:n-2]
#corr = np.corrcoef(df_cor.T)
upp_mat = np.triu(df_cor.corr())
sns.heatmap(df_cor.corr(), vmin = -1, vmax = +1, annot = True, cmap = 'coolwarm', mask = upp_mat)
<Axes: >
dftotal.iloc[:,1:14].describe()
| GSPC.Adjusted | AAPL.Adjusted | AMZN.Adjusted | MSFT.Adjusted | TSLA.Adjusted | GOOG.Adjusted | GOOGL.Adjusted | NVDA.Adjusted | BRK.B.Adjusted | META.Adjusted | UNH.Adjusted | JNJ.Adjusted | PG.Adjusted | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| count | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 | 2667.000000 |
| mean | 2685.892010 | 59.416211 | 70.665662 | 113.777254 | 70.480455 | 57.459379 | 57.640025 | 62.870372 | 185.321669 | 148.262557 | 212.637677 | 110.799351 | 88.421621 |
| std | 912.201173 | 49.474965 | 53.398954 | 91.487247 | 100.416735 | 35.656729 | 35.156219 | 74.184491 | 65.377359 | 86.401787 | 144.348470 | 35.494257 | 31.780452 |
| min | 1278.040039 | 12.046193 | 10.411000 | 21.689646 | 1.740000 | 13.924059 | 13.990240 | 2.610562 | 78.830002 | 17.730000 | 43.228558 | 45.895130 | 42.930496 |
| 25% | 1996.284973 | 22.743071 | 19.154250 | 39.305872 | 13.896667 | 28.364902 | 28.747622 | 4.884724 | 135.935005 | 78.385002 | 96.489235 | 81.087364 | 63.413342 |
| 50% | 2496.840088 | 37.335831 | 50.503502 | 71.427315 | 18.148666 | 48.944500 | 49.808498 | 37.803696 | 178.570007 | 143.800003 | 183.339828 | 112.401558 | 73.754150 |
| 75% | 3257.575074 | 86.821614 | 100.522503 | 191.692101 | 66.203666 | 73.241749 | 73.216999 | 92.120262 | 218.050003 | 190.404998 | 287.365905 | 136.825287 | 115.653481 |
| max | 4796.560059 | 180.683853 | 186.570496 | 339.075562 | 409.970001 | 150.709000 | 149.838501 | 333.350769 | 359.570007 | 382.179993 | 551.479736 | 181.108810 | 159.209808 |
El mapa de correlación nos muestra que al incluir los precios de las acciones de empresas como Microsoft (MSFT) y Berkshire Hathaway (BRK.B) en el modelo para predecir el S&P500, podríamos obtener una mejor precisión pues dado que la correlación indica una relación lineal entre dos variables, los precios de las acciones de MSFT, BRK.B y el S&P500. Como vemos, encontrmaos una correlación positiva cercana a 1 lo que puede interpretarse como que los precios de las acciones se mueven en la misma dirección que el S&P500.
De otro lado, al revisar el portafolio de acciones, observamos que el precio promedio de estas dos acciones está dentro de las acciones con los mayores precios de las seleccionadas, lo cual, vinculándolo a la forma como se calcula el S&P500 hace que estas acciones tengan una alta improtancia en el cálculo del resultado final del S&P. Si bien, hay otras acciones con mayores precios, no tienen una correlación tan alta con el índice razón por la cual no se incluirán.
Dado que el S&P500 es un índice ponderado de capitalización bursátil en donde se encuentran las 500 empresas de EEUU con gran capitalización, las empresas con mayor capitalización bursátil deben tener una mayor influencia en el resultado del índice. Adicionalmente, son empresas que pueden ser potencialmente líderes en el mercado y bien posicionadas, lo cual puede marcar tendencias en el mercado.
Por tanto, vamos a explorar diferentes posibilidadades para correr el modelo con menos variables exógenas:
dftrain_filtrado = dftrainall_est[['GSPC.ReturnSuavizado','MSTF.ReturnSuavizado','BRK.B.ReturnSuavizado','VIX.ReturnSuavizado']]
dftest_filtrado = dftestall[['GSPC.ReturnSuavizado','MSTF.ReturnSuavizado','BRK.B.ReturnSuavizado','VIX.ReturnSuavizado']]
# Configuramos la función de pronóstico de forecaster con la mejor combinación obteida de hiperparámetros.
# ==============================================================================
### Aquí debe ir tu codigo,
regressor = RandomForestRegressor(max_depth= 10, n_estimators=500 , random_state=123)
### Hasta aquí modificas
forecaster2 = ForecasterAutoreg(
regressor = regressor,
lags = 5
)
forecaster2.fit(y = dftrain_filtrado.iloc[:,0], exog = dftrain_filtrado.iloc[:,1:])
# Pronóstico, en este caso elegimos una ventana de 20 días.
# ==============================================================================
pred_filtrado = forecaster2.predict_interval(steps=steps, exog=dftest_filtrado.iloc[:,1:],interval = [1, 99],n_boot = 500)
grafico_intervalos(dftest_filtrado.iloc[:,0],pred_filtrado)
desempeño(dftest_filtrado,pred_filtrado)
EL MSE del modelo en test es de: 0.14504318029190116 EL RMSE del modelo en test es de: 0.3808453495736835 EL R2 del modelo en test es de: 0.8906032684093637 EL MAE del modelo en test es de: 0.2971541961693256
Al evaluar las diferencias en las métricas encontramos que el R2 y el mse no cambian significativamente con el ajuste de las variables, por lo que continuaremos con la experimentación para reducir al máximo la cantidad de variables a utilizar.
En el siguiente experimento, quitaremos el VIX. Al ser el VIX un índice de volatilidad, en tiempo real, que se emplea como referente para cuantificar las expectativas del mercado en relación con la volatilidad implícita del S&P 500 (SPX) durante los siguientes 30 días, tener sus valores futuros implica realizar un proceso similar al que estamos haciendo para predecir el S&P500, por lo cual analizaremos el comportamiento de este modelo sin este índice para evaluar sus resultados:
dftrain_sinvix = dftrain[['GSPC.ReturnSuavizado','MSTF.ReturnSuavizado','BRK.B.ReturnSuavizado']]
dftest_sinvix = dftest[['GSPC.ReturnSuavizado','MSTF.ReturnSuavizado','BRK.B.ReturnSuavizado']]
dftestall_sinvix = dftestall[['GSPC.ReturnSuavizado','MSTF.ReturnSuavizado','BRK.B.ReturnSuavizado']]
# Configuramos la función de pronóstico de forecaster con la mejor combinación obteida de hiperparámetros.
# ==============================================================================
### Aquí debe ir tu codigo,
regressor = RandomForestRegressor(max_depth= 10, n_estimators=500 , random_state=123)
### Hasta aquí modificas
forecaster_sinvix = ForecasterAutoreg(
regressor = regressor,
lags = 5
)
forecaster_sinvix.fit(y = dftrain_sinvix.iloc[:,0], exog = dftrain_sinvix.iloc[:,1:])
# Pronóstico, en este caso elegimos una ventana de 20 días.
# ==============================================================================
pred_sinvix = forecaster_sinvix.predict_interval(steps=steps, exog=dftest_sinvix.iloc[:,1:],interval = [1, 99],n_boot = 500)
grafico_intervalos(dftest_sinvix.iloc[:,0],pred_sinvix)
desempeño(dftest_sinvix,pred_sinvix)
EL MSE del modelo en test es de: 0.22132873517154628 EL RMSE del modelo en test es de: 0.47045588015407597 EL R2 del modelo en test es de: 0.8342890621228852 EL MAE del modelo en test es de: 0.3631235752483951
Aunque si perdemos desempeño, encontramos que la predicción del VIX resulta tan compleja como la predicción del S&P500, razón por la cual se contempla omitir esta variable.
Sin embargo, antes de tomar la decisión veamos qué tanto podría aportar al modelo únicamente el VIX:
dftrain_filtrado2 = dftrainall_est[['GSPC.ReturnSuavizado','VIX.ReturnSuavizado']]
dftest_filtrado2 = dftestall[['GSPC.ReturnSuavizado','VIX.ReturnSuavizado']]
# Configuramos la función de pronóstico de forecaster con la mejor combinación obteida de hiperparámetros.
# ==============================================================================
### Aquí debe ir tu codigo,
regressor = RandomForestRegressor(max_depth= 10, n_estimators=500 , random_state=123)
### Hasta aquí modificas
forecaster_vix = ForecasterAutoreg(
regressor = regressor,
lags = 5
)
forecaster_vix.fit(y = dftrain_filtrado2.iloc[:,0], exog = dftrain_filtrado2.iloc[:,1:])
# Pronóstico, en este caso elegimos una ventana de 20 días.
# ==============================================================================
pred_filtrado2 = forecaster_vix.predict_interval(steps=steps, exog=dftest_filtrado2.iloc[:,1:],interval = [1, 99],n_boot = 500)
grafico_intervalos(dftest_filtrado2.iloc[:,0],pred_filtrado2)
desempeño(dftest_filtrado2,pred_filtrado2)
EL MSE del modelo en test es de: 0.7542956626596483 EL RMSE del modelo en test es de: 0.868501964683816 EL R2 del modelo en test es de: 0.4310833506139946 EL MAE del modelo en test es de: 0.703442469348462
Este resultado nos muestra que si bien el VIX es importante, nos arroja una mayor precisión agregar los predictores de las acciones de MSTF y BRK.B, por lo cual la decisión final de los analistas será trabajar en dos modelos adicionales para predecir el comportamiento de estas dos acciones y omitir los datos del VIX como variable de entrada para el modelo
Corramos el modelo con el 100% de los datos de entrenamiento ahora:
dftrainall_sinvix = pd.concat([dftrain_sinvix,dftest_sinvix])
dftestall_sinvix = dftestall[['GSPC.ReturnSuavizado','MSTF.ReturnSuavizado','BRK.B.ReturnSuavizado']]
# Configuramos la función de pronóstico de forecaster con la mejor combinación obteida de hiperparámetros.
# ==============================================================================
### Aquí debe ir tu codigo,
regressor = RandomForestRegressor(max_depth= 10, n_estimators=500 , random_state=123)
### Hasta aquí modificas
forecaster_sinvix = ForecasterAutoreg(
regressor = regressor,
lags = 5
)
forecaster_sinvix.fit(y = dftrainall_sinvix.iloc[:,0], exog = dftrainall_sinvix.iloc[:,1:])
# Pronóstico, en este caso elegimos una ventana de 20 días.
# ==============================================================================
pred_sinvixfinal = forecaster_sinvix.predict_interval(steps=steps, exog=dftestall_sinvix.iloc[:,1:],interval = [1, 99],n_boot = 500)
grafico_intervalos(dftestall_sinvix.iloc[:,0],pred_sinvixfinal)
desempeño(dftestall_sinvix,pred_sinvixfinal)
EL MSE del modelo en test es de: 0.18220597162151003 EL RMSE del modelo en test es de: 0.4268559143569526 EL R2 del modelo en test es de: 0.862573767814698 EL MAE del modelo en test es de: 0.37095388952853675
Tomemos entonces las predicciones obtenidas con el modelo sin el VIX para encontrar los resultados del modelo en retornos, recordemos que escalamos los datos y debemos volverlos a sus dimensiones originales:
# Construimos una matriz igual a la del conjunto de test para cada uno de los 3 intervalos de predicción
# ==============================================================================
#Matriz para predicciones
pred_filtrado_return = pd.DataFrame(0, index=pred_sinvixfinal.index, columns=dftestall.columns)
pred_filtrado_return['GSPC.ReturnSuavizado'] = pred_sinvixfinal['pred']
#Matriz para intervalo inferior
predlow_filtrado_return = pd.DataFrame(0, index=pred_sinvixfinal.index, columns=dftestall.columns)
predlow_filtrado_return['GSPC.ReturnSuavizado'] = pred_sinvixfinal['lower_bound']
#Matriz para intervalo superior
predup_filtrado_return = pd.DataFrame(0, index=pred_sinvixfinal.index, columns=dftestall.columns)
predup_filtrado_return['GSPC.ReturnSuavizado'] = pred_sinvixfinal['upper_bound']
#Transformamos los datos a la representación original.
pred_returns1 = pd.DataFrame(scaler.inverse_transform(pred_filtrado_return, copy=None))[0]
predlow_returns1 = pd.DataFrame(scaler.inverse_transform(predlow_filtrado_return, copy=None))[0]
predup_returns1 = pd.DataFrame(scaler.inverse_transform(predup_filtrado_return, copy=None))[0]
#Construimos un dataframe con las predicciones en la representación original
pred_final_1 = pd.concat([pred_returns1,predlow_returns1,predup_returns1],axis=1)
pred_final_1.columns = ["pred","lower_bound","upper_bound"]
pred_final_1.index = pred_sinvixfinal.index
# Importamos la información sin estandarizar de los retornos que tenemos en la predicción
# ==============================================================================
return_real = pd.DataFrame(df2.iloc[2792:2812,0])
return_real.head()
| GSPC.ReturnSuavizado | |
|---|---|
| Date | |
| 2023-02-01 | 1.039806 |
| 2023-02-02 | 1.459238 |
| 2023-02-03 | -1.040859 |
| 2023-02-06 | -0.615939 |
| 2023-02-07 | 1.279036 |
# Calculamos las métricas
# ==============================================================================
grafico_intervalos(return_real.iloc[:,0],pred_final_1)
desempeño(return_real,pred_final_1)
EL MSE del modelo en test es de: 0.14076582767600526 EL RMSE del modelo en test es de: 0.3751877232479832 EL R2 del modelo en test es de: 0.8625737678146981 EL MAE del modelo en test es de: 0.32605228265810177
dftestall['GSPC.ReturnSuavizado'].max(),dftestall['GSPC.ReturnSuavizado'].min()
(1.5971519617016467, -2.3663176740723677)
Finalmente, vamos a entrenar el modelo con todos los datos disponibles para futuras predicciones:
df_completo = pd.concat([dftrainall_sinvix,dftestall_sinvix])
# Configuramos la función de pronóstico de forecaster con la mejor combinación obteida de hiperparámetros.
# ==============================================================================
### Aquí debe ir tu codigo,
regressor = RandomForestRegressor(max_depth= 10, n_estimators=500 , random_state=123)
### Hasta aquí modificas
forecaster_final = ForecasterAutoreg(
regressor = regressor,
lags = 5
)
forecaster_final.fit(y = df_completo.iloc[:,0], exog = df_completo.iloc[:,1:])
import joblib
joblib.dump(scaler, 'scaler.gz')
['scaler.gz']
# Guardar modelo
#=====================================================================
save_forecaster(forecaster_final, file_name='forecaster.py', verbose=False)
load_forecaster('forecaster.py')
=================
ForecasterAutoreg
=================
Regressor: RandomForestRegressor(max_depth=10, n_estimators=500, random_state=123)
Lags: [1 2 3 4 5]
Transformer for y: None
Transformer for exog: None
Window size: 5
Weight function included: False
Exogenous included: True
Type of exogenous variable: <class 'pandas.core.frame.DataFrame'>
Exogenous variables names: ['MSTF.ReturnSuavizado', 'BRK.B.ReturnSuavizado']
Training range: [Timestamp('2012-05-21 00:00:00'), Timestamp('2023-02-28 00:00:00')]
Training index type: DatetimeIndex
Training index frequency: B
Regressor parameters: {'bootstrap': True, 'ccp_alpha': 0.0, 'criterion': 'squared_error', 'max_depth': 10, 'max_features': 1.0, 'max_leaf_nodes': None, 'max_samples': None, 'min_impurity_decrease': 0.0, 'min_samples_leaf': 1, 'min_samples_split': 2, 'min_weight_fraction_leaf': 0.0, 'n_estimators': 500, 'n_jobs': None, 'oob_score': False, 'random_state': 123, 'verbose': 0, 'warm_start': False}
Creation date: 2023-05-07 16:26:00
Last fit date: 2023-05-07 16:26:09
Skforecast version: 0.7.0
Python version: 3.9.13
Forecaster id: None
=================
ForecasterAutoreg
=================
Regressor: RandomForestRegressor(max_depth=10, n_estimators=500, random_state=123)
Lags: [1 2 3 4 5]
Transformer for y: None
Transformer for exog: None
Window size: 5
Weight function included: False
Exogenous included: True
Type of exogenous variable: <class 'pandas.core.frame.DataFrame'>
Exogenous variables names: ['MSTF.ReturnSuavizado', 'BRK.B.ReturnSuavizado']
Training range: [Timestamp('2012-05-21 00:00:00'), Timestamp('2023-02-28 00:00:00')]
Training index type: DatetimeIndex
Training index frequency: B
Regressor parameters: {'bootstrap': True, 'ccp_alpha': 0.0, 'criterion': 'squared_error', 'max_depth': 10, 'max_features': 1.0, 'max_leaf_nodes': None, 'max_samples': None, 'min_impurity_decrease': 0.0, 'min_samples_leaf': 1, 'min_samples_split': 2, 'min_weight_fraction_leaf': 0.0, 'n_estimators': 500, 'n_jobs': None, 'oob_score': False, 'random_state': 123, 'verbose': 0, 'warm_start': False}
Creation date: 2023-05-07 16:26:00
Last fit date: 2023-05-07 16:26:09
Skforecast version: 0.7.0
Python version: 3.9.13
Forecaster id: None